

adjust_pval_step_down_fun = function(p_val_original, p_val_RI_mat){
  # p_val_original: unadjusted p-values using entire sample
  # p_val_RI_mat: matrix where each row is one vector of p-vals from the RI procedure
  
  # This auxiliary data frame allows me to sort the RI pvalues based on ascending order of the original p-values
  df.aux = tibble(p_val_original, pos = seq(1, length(p_val_original))) %>% 
    arrange(p_val_original)
  
  p_vals_original_s = df.aux[["p_val_original"]]
  pos.vec = df.aux[["pos"]]
  
  # Sort RI matrix so it follows the order of the 
  p_val_RI_mat_s = p_val_RI_mat[, pos.vec]
  
  ## Step down procedure
  
  # Matrix that captures for each iteration of the RI procedure whether we reject each of the outcomes
  
  R.iter.mat = matrix(NA_real_, nrow = nrow(p_val_RI_mat_s), ncol = ncol(p_val_RI_mat_s))
  for(iter in 1:nrow(p_val_RI_mat_s)){
    
    p.ri.iter = p_val_RI_mat_s[iter,]
    # vector to store RI p-values for each iteration of RI procedure after imposing monotonicity
    p.ri.mon.iter = rep(NA_real_, ncol(p_val_RI_mat_s))
    p.ri.mon.iter[1] = min(p.ri.iter)
    for(r in 2:ncol(p_val_RI_mat_s)){
      p.ri.mon.iter[r] = min(p.ri.iter[-seq(1,r-1)])
    }
    
    # Compute whether p_r** < p_r (where )
    R.iter.mat[iter, ] = (p.ri.mon.iter < p_vals_original_s)
  }
  # Number of times we reject for each of the outcomes
  S.r = colSums(R.iter.mat)
  
  # Step 6: Calculate p-value without imposing monotonicity 
  p_fwer = S.r/nrow(p_val_RI_mat_s)
  
  # Step 7: Impose monotonicity in order to p_fwer
  p_fwer_mon = p_fwer
  p_fwer_mon[1] = min(p_fwer)
  for(r in 2:ncol(p_val_RI_mat_s)){
    p_fwer_mon[r] = min(p_fwer[-seq(1,r-1)])
  }
  
  # sort p_fwer_mon so it matches the order the original pvalue was inputted
  df.aux = tibble(
    p_fwer_mon = p_fwer_mon,
    pos = pos.vec
  ) %>% 
    arrange(pos)
  
  p_fwer_mon_ori_order = df.aux$p_fwer_mon
  
  return(p_fwer_mon_ori_order)
  
}
